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Abstract 

We consider the persistence and stability of dark solitons in the Gross-Pitaevskii (GP) equation 
with a small decaying potential. We show that families of black solitons with zero speed originate 
from extremal points of an appropriately defined effective potential and persist for sufficiently small 
strength of the potential. We prove that families at the maximum points are generally unstable 
with exactly one real positive eigenvalue, while families at the minimum points are generally un- 
stable with exactly two complex-conjugated eigenvalues with positive real part. This mechanism 
of destabilization of the black soliton is confirmed in numerical approximations of eigenvalues of 
the linearized GP equation and full numerical simulations of the nonlinear GP equation with cubic 
nonlinearity. We illustrate the monotonic instability associated with the real eigenvalues and the 
oscillatory instability associated with the complex eigenvalues and compare the numerical results of 
evolution of a dark soliton with the predictions of Newton's particle law for its position. 



1 Introduction 

Dark solitons are solutions of nonlinear PDEs in the space of one dimension with non-zero boundary 
conditions and non-zero phase shift. They are represented by a family of traveling wave solutions 
extending from the limit of zero speed (so-called black solitons) to the limit of sound speed (so-called 
grey solitons). Prom a physical point of view, dark solitons are waves in defocusing nonlinear systems 
which move along the modulationally stable continuous-wave background. 

The original interest in studies of dark solitons emerged, roughly, two decades ago in the context 
of nonlinear optics, where dark solitons provide modulations of the light intensity of an optical beam 
traveling in a planar waveguide (22]. The main model for dark solitons in nonlinear optics is the 
generalized nonlinear Schrodinger (NLS) equation 

iu t = ~^u xx + f(\u\ 2 )u, (1.1) 

where u(x,t) : R x R + C and f{q) : R + — * R. We assume that f{q) is a smooth monotonically 
increasing function on q £ I C R + . Particular examples of f(q) include the cubic NLS with f(q) = q, 
the cubic-quintic NLS with f(q) = aq + q 2 , a S R and the saturable NLS with / = —1/(1 + f3q), 
P G R+. 
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Among others, several analytical results were important in the development of dark solitons in recent 
years: perturbation theory based on renormalized power |24| and momentum orbital stability of 
dark solitons ^ |2] , completeness of eigenfunctions in the cubic NLS [SJ I15j , inverse scattering for the 
vector cubic NLS equation construction of the Evans function for dark solitons in the perturbed 
cubic NLS |19j , asymptotic analysis of the radiation and dynamics of dark solitons |261 1271 I35| , and 
spectral analysis of transverse instabilities of one-dimensional dark solitons |25l I36| . 

Subsequently, rigorous analysis of the existence and stability of dark solitons was developed in the last 
decade based on the earlier physical literature. In particular, Zhidkov |46j proved local existence of the 
Cauchy problem and stability of kink solutions, de Bouard [IJ proved spectral and nonlinear instability 
of stationary bubbles (black solitons with zero phase shift), Lin |2H] proved the criterion for orbital 
stability and instability of dark solitons (for non-zero velocities), Maris studied bifurcation of dark 
solitons (for non-zero velocities) in the presence of delta-function potential and its generalizations, and 
Di Menza and Gallo [5j proved recently the stability criterion for kinks (black solitons with a non-zero 
phase shift). Extensions of the existence and stability theory to the generalized NLS equation in two 
and higher dimensions were also developed by J.C. Saut and his co-workers. 

While many mathematical results are now available for solutions of the generalized NLS equation 
Ql.l|) with non-zero boundary conditions and non-zero phase shift, dark solitons have suffered a decreas- 
ing popularity in the context of nonlinear optics. This is not only because they possess infinite energy 
due to non-zero boundary conditions but also because it is difficult from the experimental point of view 
to separate the effects of the dark soliton dynamics and the background dynamics. Also no engineering 
or commercial applications of dark solitons have been proposed, to the best of our knowledge. 

Nevertheless, the interest in these nonlinear waveforms has recently been rejuvenated by the rapid 
development of a new area of physics, namely the field of Bose-Einstein condensates (BECs) |39|I4()|. In 
the latter setting, dark solitons typically move along the localized ground state trapped by the external 
potentials [12]. The main model for BECs is a modification of the NLS equation with an external 
potential, which is called the Gross-Pitaevskii (GP) equation, 

iu t = -~u xx + f(\u\ 2 )u + eV(x)u, (1.2) 

where e E R is the strength of the potential V{x) and V(x) : R — ► R is assumed to be a smooth function 
satisfying one of the three properties: 

(i) V(x) is bounded and decaying, e.g. 

3C>0, k>0: \V(x)\ < Ce~ K \ x \ VxEl (1.3) 

(ii) V(x) is bounded but non-decaying, e.g. V(x + d) = V(x), i£K with period d > 

(iii) V(x) is unbounded, e.g. V(x) = x 2 + V(x), where V(x) is bounded on x 6 R. 

The last two cases are of particular interest in the context of Bose-Einstein condensates, where dynamics 
of dark solitons (localized dips in the nonlinear ground state trapped by the potential V(x)) was studied 
in many recent papers, see 0133] for surveys of results. Although cases (ii) and (iii) have been initially 
our primary motivations, this paper covers only the case (i) when V(x) is bounded and exponentially 
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decaying as in (|1.3|) . In particular, we consider the class of symmetric potentials V(—x) = V(x) with 
two examples 

Vi(x) = -sech 2 , V 2 (x) = xV kN , x£R. (1.4) 

The other cases (ii) and (iii) of the potential V(x) will be the subject of forthcoming studies. 

We note in passing that while the potentials (|1.4|) are, perhaps, less customary than the standard 
magnetic (parabolic) and optical lattice (periodic) potentials |21j . they are nonetheless still quite phys- 
ically relevant. In particular, the potential V\{x) corresponds to a red-detuned laser beam potential, 
analogous to the one used in [HSUE]. The potential V2(x) represents an all-optically trapped BEC, as 
modeled in [7j and experimentally implemented in [3]. 

The strategy of our work is to exploit solutions of the GP equation (|1.2j) in the limit of small 
strength e. Starting with the limit e = 0, where both existence and stability of dark solitons with 
non-zero and zero speeds are known from the analysis of the NLS equation (|1.1|) . we shall use the 
method of Lyapunov-Schmidt reductions and obtain existence and stability results for small e. From 
the technical point of view, we use local bifurcation analysis of solutions of the ODEs with non-zero 
boundary conditions similarly to [SO], persistence analysis of isolated eigenvalues in the problems with 
small potentials similarly to 18 , the count of eigenvalues in the generalized eigenvalue problem for 
self-adjoint operators with no spectral gaps similarly to 0, and the construction of L 2 eigenfunctions 
of the stability problem with fast and slow exponential decay similarly to |SH]- Since our starting point 
is the case of e = 0, we will also give alternative rigorous proofs to the existence and stability of black 
solitons in the NLS equation Ijl.lJ) . which complement the recent work of |31j . 

Our main results are listed as follows. 



(i) Let u = 4> (x - s ) e - iujt+w be a solution of the NLS equation JlTJ with s G R, 9 G R, u = f(q ), 
qo G T C M+ and (fto( x ) -R ^R such that 4>o(x) — » ±y/q~o as x — > ±oo exponentially fast. Let sq 
be the root of the function M'(s), such that 

M'(so) = / V'{x) [go - (f> 2 {x - s )] dx = (1.5) 

JR. 

and 

M"(s ) = / V"(x) [q - ct>l(x - s )] dx + 0. (1.6) 
JR 

Then, there exists a unique continuation of the solution u = (p t (x — s e )e~ tujt+tS of the GP equation 
(jl.2j) with V(x) in (|1.3|) and e sufficiently small, where 4> e : R i— > R and 4> e (x) — > ±^/q~Q as x — > ±oo 
exponentially fast, such that (fi t (x) and s e are e-close to 4>o(x) and sq in L°°-norm. 

(ii) Let the solution (f>o(x) be spectrally stable in the time evolution of the NLS equation (|1.1|) . 
Then, the solution (p e {x) is spectrally unstable in the time evolution of the GP equation (|1.2j) 
for sufficiently small e with exactly one real positive eigenvalue if M"(so) < and exactly two 
complex-conjugate eigenvalues with positive real part if M"(sq) > 0. 

(iii) Let the function uo(x) : R i—* R be close to the solution 4> e (x) in L°°-norm and satisfy the same 
boundary conditions uq(x) — > ±^/qo as x — > ±oo. Then, the solution u(x, t) : R X i— > C of the 
Cauchy problem to the GP equation (jl.2fl with it(x,0) = uq(x) remains close to the dark soliton 
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4> € {x — s(t)) with varying coordinate s(t) for < t < T such that s(t) solves the Newton's particle 
equation 

Ho's - e\ M"(s)s = -eM'(s), < t < T, (1.7) 

where M(s) is the effective potential in (|1.5|) and (jl.6j) . (fio,Xo) are constants representing the 
soliton's mass and gain terms, and the time T > is of the order of 0(l/e). 

Statement (i) is formulated and proved in Section 2 (see Theorem l2.12j) . Statement (ii) is formulated 
and proved in Section 3 (see Theorems 13.111 13.141 and I3.15j) under some non-degeneracy assumptions 
(see Cor pilar v 13 . 1 71 and Remark l3.18|) . The two complex-conjugate eigenvalues with positive real part in 
the statement (ii) for M"(so) > follows from the linearized version of the Newton's particle equation 
Q1.7|) with fiQ > and Ao > 0, which is rigorously derived in Section 4 (see Theorem I4.11|) . Constants 
(/Uo, Ao) are defined by E,emark l4.131 The linearization around a black soliton is investigated in Section 
5 using numerical bifurcation theory and approximations of eigenvalues. Statement (iii) is a conjecture 
which is tested against appropriately crafted numerical experiments in Section 6. The summary and 
open problems are discussed in the concluding Section 7. 

2 Existence analysis of dark and black solitons 

We first consider the family of traveling solutions of the NLS equation After we formulate the 

conditions for existence of dark and black solitons, we address persistence of stationary solutions in 
the Gross-Pitaevskii equation (|1.2|) for small e and the potential V(x) in (jl.3fl . Applications of the 
persistence analysis to the two potentials (|1.4j) indicate that the families of black solitons bifurcate 
from the extremal points of the effective potential M{s) in (|1.5|) - (|1.6|) . 

Definition 2.1 A dark soliton is the traveling solution of the NLS equation in the form: 

u(x,t) = U{x- vt)e~ iuJ \ U(z) = <$>{z)e ie{z \ z = x - vt, (2.1) 

where U : R \—* C, R + , and 0:Ih [— tt,tt] are smooth functions of their arguments, which 

converge exponentially fast to the boundary conditions 

lim $(z) = y/%, lim &(z) = G±. (2.2) 

Here (u,v) G V C R 2 , q G 1 C M+ and G± G [— 7r, 7r] are parameters of the solution. Moreover, the 
functions <&(z) and Q(z) can be chosen to satisfy the normalization conditions $'(0) = and + = 0. 
Additionally, we require that &(z) < ^/qo on z G R. 

Remark 2.2 The normalization conditions for &(z) and @(z) use the gauge [u(x,t) — > u(x,t)e ld , 
\/0 G M] and translational [u(x,t) — > u(x — s,t), Vs G M] invariance of the NLS equation, while the 
linear growth of @(z) in z is excluded by the Galileo invariance [u(x, t) — > u(x — kt, t)e tkx ~' lk */ 2 , V/c G M]. 



Theorem 2.3 Let f(q) be C 1 (IR + ) and fix qo G R+ such that f'(qo) > and c = ^qof'iqo) > 0. 
The dark soliton U{z) of Definition \2.1\ exists if u) = f(qo) and v G (— c, 0) U (0, c). Moreover, for 
these solutions, &(z) has a global minimum at z = with < $(0) < ^/q~Q and Q(z) is monotonically 
decreasing for v > and increasing for v < 0. 
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Proof. It follows immediately that the function U(z) satisfies the second-order ODE: 

-ivU' + ~U" + (u-f(\U\ 2 ))U = 0, (2.3) 

while the functions &(z) and @(z) satisfy the ODE system in the hydrodynamic form: 

$"_(e')2$ + 2t ,e'$ + 2(a;-/($ 2 ))$ = (2.4) 

(0'$ 2 -v<& 2 )' = (2.5) 

Integrating the ODE (|2.5|) under the boundary conditions (|2.2|) . we obtain 

*-4^fi. (2.6) 

As a result, the ODE (|2.4j) reads as follows 

$" + 2{lo - /($ 2 ))$ + v 2 ^~ q ° = 0. (2.7) 

The equilibrium point = exists if and only if u = f(qo) and it is a non-degenerate hyperbolic 
point if and only if v 2 < c 2 = qof'(qo). Integrating the second-order ODE ()2.7jl subject to the boundary 
conditions (|2.2j) . we obtain the quadrature 

($') 2 ~ 2IV ($ 2 ) + v 2 ^ ^" 2 9o)2 = 0, (2.8) 



where 



W(q)= I" (f(q )-f(q))dq. (2.< 
J qo 

,2/ 



The non-degenerate turning point < <I>o < ^/qo of the effective potential — 2W{q) + v (q — qo)/q exists 
for any < v 2 < c 2 , such that the trajectory from the hyperbolic point turns at <I>(0) = <&o and returns 
back to the hyperbolic point forming a homoclinic orbit, thus proving the statement. □ 

Remark 2.4 It is not difficult to prove from the ODE analysis that U{z) is in fact C°°(R) for v G 
(— c, 0) U (0, c) and U (z) converges to the limiting values ^/qoe' l@± as z — > ±oo at the exponential 
rate 0(e~ v ' c ' 2 ~ v2z ). Although Lemma 12.31 formulates only the sufficient condition for existence of dark 
solitons, it is possible to consider all other solutions of the ODE system for <&(z) and @(z). The case 
v 2 > c 2 corresponds to the elliptic point = yfq$ and has no homoclinic orbits at all. The marginal 
case v 2 = c 2 may admit an algebraically decaying homoclinic orbit, which is excluded by Definition 
12.11 In the case v 2 < c 2 , other homoclinic orbits can exist with &(z) > ^/q$ but they are also excluded 
by Definition 12.11 The only case, which may be included into the necessary and sufficient condition of 
existence of dark soliton, is the limit v — > 0. This limit leads to two different kinds of solutions which 
are classified in the following definition. 



Definition 2.5 The limiting solution (fo[x) = \\mU{x) in the family of dark solitons of Definition ^. l\ 

is said to be the black soliton if <fio(x) is a real-valued smooth function on x G R. The black soliton is 
called the bubble if <po(~ x ) = (ftoix) with < (j>o(0) < Jqf) and lim (fro(x) = Jqo, while it is called the 

x— >±oo 

kink if (J)q(— x) = — <f>o(x) with 0o(O) = and lim <fto{x) = ±y/qo, where 4>o( x ) converges exponentially 

x— >±oo 

fast to the limits. 
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Theorem 2.6 Let f(q) satisfy the same conditions as in Lemma \2.cA and W{q) be defined by \2. y\) . 
The kink solution of Definition \2. 51 exists if and only if to = f(qo), v = 0, and W{q) > for all 
< q < qo- The bubble solution of Definition \2. 51 exists if and only if to = f(qo), v = 0, and there 
exists the largest value q\ G [0,go) that satisfies W{q\) = and W'(qi) > 0. 

Proof. By Lemma \'2. 31 we have < <&(z) < yf% and Q'(z) < on z G K for v > and hence no kink or 
bubble solutions may exist for v > 0. Let us consider the real- valued solution (fro(x) °f the second-order 
ODE 

^ + (w - = 0. (2.10) 

The equilibrium points (fro = ±^/qo exist if and only if lo = f(qo) and they are non-degenerate hyperbolic 
points. Integrating the second-order ODE (|2.10|) under the boundary conditions in Definition 12.51 we 
obtain 

(<// ) 2 - 2W{4>1) = o, 

where W(q) is defined by (|2.9|) , No turning point of W(q) exists if W(q) > for any q G [0, go); and 
the outgoing trajectory from the hyperbolic point (fro = ^/qo connects the incoming trajectory to the 
hyperbolic point (fro = —^/q~o an d vice versa forming a pair of heteroclinic orbits. If there exists a non- 
degenerate largest turning point q\ of W(q) on q\ G [0, qo), such that W(q\) = and W'(q\) > 0, then 
the trajectory from the hyperbolic point 4>q = ^/qo turns at (po = yj~q~\ and returns back to the point 
00 = \/Qo forming a homoclinic orbit. If the point q\ is degenerate, such that W^f/i) = V^'((7i) = 0, 
there exists a front solution <po{x) from lim (fio( x ) = \/Qi to lim ^o(^) = Jqo and vice versa forming 

a pair of fronts, which are excluded by Definition 12.51 

It remains to prove that the family of dark solitons U(x) for < v < c converges to the black soliton 
(fro (x) as v I 0. The proof follows from the quadrature ()2.8j) . There exist unique classical solutions (fr±{z) 
on z G R± for any v > with $±(0) = ^fql found from the largest root of — 2W{q*) +v 2 {q* — qo)/q* = 
on q G [0,(/o). In the limit v [ 0, the root q* converges to if W(0) > and no other roots of W(q) 
exists on q G [0, go). Otherwise, the root converges to 51, where q\ is the largest root of W(q) on 
q G [0,qo)- It follows from the ODE (|2.6|) in the limit v j that Q(z) is piecewise constant function on 
z£l with a possible jump discontinuity at z = 0. 

Let us write x = z for u = 0. In the case g* [ as i> j 0, the two smooth solutions &±(x) are 
glued into one smooth real- valued solution 4>o( x ) if and only if Q(x) = on x G M+ and O(x) = 7r on 
x G R_ (under the normalization + = 0). This limiting solution becomes the kink of Definition 12.51 
In the case q* [ q± as v [ and q± is a non-degenerate root of W(q), the two smooth solutions 3>±(x) 
are glued into a smooth real- valued solution <fro(x) if and only if O(x) = on x G M (under the same 
normalization). This limiting solution becomes the bubble of Definition 12.51 □ 

Remark 2.7 Due to their potential stability in the time evolution of the NLS equation only 
kinks are considered in the GP equation (|1.2|) for sufficiently small e. The bubbles are always unstable 
in the time evolution of the NLS equation (jl.ll) 4 . 

Definition 2.8 The kink mode is the stationary solution of the GP equation \1.2\) in the form: 

u(x, t) = <t> e (x)e- if ^ t+ie , flel, 

where 4> t (x) is a real-valued smooth function on x G M., which converges exponentially fast to the 
boundary conditions lim (fr t (x) = ±y/qo, for any e G M. 
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Lemma 2.9 Let 4> t (x) be the kink mode of Definition \2.8\ and V{x) be the potential satisfying hl.'J\) . 
Then, for any e ^ 0, 

r V{x) [q - 4> 2 e (x)] dx = 0. (2.11) 



Proof. The stationary solutions of Definition 12 . 81 satisfy the second-order ODE: 

\€ + (/(go) " f(<%)) & = eV(x)tj> e , (2.12) 
which is generated by the Hamiltonian function 

E(<f> e ,<t>' £ ,x) = i(0' £ ) 2 - W{<%) + € -V{x) [q - , 

where W(q) is given by Q2.9B . Therefore, the change of E{(f> e {x), (j>' £ (x),x) at the classical solution 4> e (x) 
of the second-order ODE (|2.12j) is given by 

dE 

— = eV'(x) [go - <t>l{x)} . 

Integrating this equation on x £ R and using the boundary conditions for <p t {x) of Definition 12.81 such 
that lim E(dJx). d[(x). x) = 0. we derive the condition (|2.11|) . □ 

x— »±oo 

Remark 2.10 If V(—x) = V{x) and 4> e (—x) = —<j) e (x) on x E R, the necessary condition \2.11\) is 
always satisfied. The center of the kink is located at x = 0, which is the minimal point of V(x) if 
V"(0) > and maximal point ifV"(0) < 0. 

Remark 2.11 The necessary condition ()2.11j) specifies restrictions on the shape of the kink mode 
4>e( x ) but does not give us any information about its existence. By using the smallness of the external 
parameter e in front of V(x), we will show that this condition is equivalent to the bifurcation equation 
in the Lyapunov-Schmidt reduction technique. A similar result for bright solitons was obtained in ^H] . 
A different role of the condition l|2.11|) was exploited in |SZj in the context of local bifurcations of small 
gap solitons in finite periodic potentials V(x). 

Theorem 2.12 Let <f>o(x) be the kink of Definition \2. 5l Let sq be a simple root of the function 

M'(s) = / V'(x) [g - <t>l{x - s)] dx, s£l, (2.13) 

such that M'(s ) = and M"(s ) / 0. Let f(g) be C 1 (R + ) and V{x) be C 2 (R) satisfying JO) . 
Then, there exists a unique continuation of (po(x — sq) to the kink mode 4> e {x — s € ) of Definition ^. 81 for 
sufficiently small e, such that 4> e (x) and s e are e-close to 4>q{x) and sq in the LP -norm. 

Proof. We use the decomposition 4> € (x) = 4>o(x — s) + cp(x,e,s) and rewrite the second-order ODE 
(|2.12|) for (f> t (x) as the root of the nonlinear operator- valued function 

F{ip, e, s) = L +i p + N(<p, s, e) + eV(x) [0 o (x -s) + <p], (2.14) 
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where L + : H 2 {W) ^ L 2 (IR) is the self- adjoint operator parameterized by s 

L+ = ~\dl + ml) - /(5b) + 2^/'(^), 
and N(ip, e, s) : H 2 (M.) <— > i? 2 (IR) is the nonlinear vector field 

iV = [/((0o + ^) 2 ) - /($) - 2^/'($)] + ^ [/(fa> + V?) 2 ) - /(0o)] . 

such that N(ip,e,s) = o(\\ip\\ H 2) as |M|#2 — > (since / S C^M^-)). Because 000*0 converges to 
iy/go as i -> ±oo exponentially fast, the essential spectrum of L + is bounded from below by 2c 2 > 
0. The operator may have isolated positive eigenvalues and no negative eigenvalues since the kernel 
L + (f)'(x — s) = is a positive definite ground state. Therefore, H 2 (K) = Ker(L + ) © Ker(L + )- L and 
the method of Lyapunov-Schmidt reductions can be applied. Projection of F(ip,e,s) onto Ker(L + )^ 
defines a unique smooth map (x,e,s) i— ► (p E Ker(L + )- L C H 2 (R) such that ||</?||#2 = 0(e) as e — > 0. 
By the Sobolev Embedded Theorem, H^IU 00 = 0(e), such that <p(x) can be decomposed as follows: 

(p = e(pi(x) + <p(x,e,s), 

where fi(x) is specified below and = o(e). Projection of F(ip,e,s) onto Ker(L + ) defines the 

bifurcation equation: 

G(e, s) = e («4 F(x)(0 o + <p)) + ($,, e, a)) = |m'(s) + G(e, s), (2.15) 

where (•,•) is the standard inner product in L 2 (M) and G(e,s) = o(e) as e — > 0. By the Implicit 
Function Theorem for the root of (|2.15j) . we obtain that s = so + s(e), where M'(sq) = 0, M"(so) ^ 0, 
and s = o(l) as e — > 0. By the Lyapunov-Schmidt Reduction Theorem for the root of (|2.14|) . a unique 
continuation of (f>o(x — sq) into (j) e (x — s t ) exists. In particular, the correction term <p\(x) satisfies the 
inhomogeneous problem 

L+<p 1 = -V(x)<l> (x-*o), ( 2 - 16 ) 

which has a unique solution cpi G Ker(L + )^ C -ff 2 (K) by the Fredholm Alternative (since M'(so) = 0). 
□ 

Corollary 2.13 (i) Let /(g) be C 2 (R + ). Then, N(<p,e,s) = 0(\\(p\\ H2 ) and G(e,s) = O(e) as e -» 0, 
suc/i £/taf s = O(e) and \\<p\\l°° = 0(e 2 ) as e ^ 0. 

(ii) Let /(g) be C°°(M + ). Then, s(e) and <p(x,e,s) are C°° -functions of e near e = 0. 

Remark 2.14 The renormalization of (p 2 t— > (^ 2 — go) is n °t needed if the potential V(a;) satisfies the 
condition 1)1. 3j) . We expect that the same quantity M(s) with the renormalization above can be useful 
to treat the other cases (ii) and (iii) of the potential term V(x). However, the method of Lyapunov- 
Schmidt reductions does not work for these cases since the essential spectrum of L + is structurally 
deformed if e 7^ and V(x) is not decaying. In this case, the root finding problem (|2.14|) is not defined 
in H 2 (M) because the term V(x)4> e (x) does not decay to zero exponentially fast as \x\ — > 00. 

Example 2.15 When the cubic NLS is considered with /(s) = s, the value go £ ^+ can be normalized 
by go = 1- in this case, the second-order ODE (|2.3j) admits an exact solution 

U(z) = kt&n)\{kz) +iv, k = \/\-v 2 , (2.17) 
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where v 6 (— 1, 1). The other components <&(z) and Q(z) of the polar form (|2,1|) can also be evaluated 
explicitly. The black soliton corresponds to the kink 4>o(x) = tanhx. When the potential V(x) is even 
V(—x) = V(x), the function M'(s) in Q2.13|) can be split into two parts: 




V'{x) [<f%(x + s) - <ffi(x - s)) dx = L'(s) - L'(-s) 



where 



L 



(s) = / V(x) [q -4>l(x- s)} dx 



(2.18) 



If V(x) is C 2 (R) and satisfies the decay condition (|1.3|) . the function L(s) is C 2 (R) and L(s) — > 
exponentially fast as |s| — ► oo. Therefore, M'(0) = -£'(0) — L'(0) = and one family of kink modes 
bifurcates from sq = 0. Additional families of kink modes of the GP equation H1.2j) may bifurcate if 
L(0) and L"(0) are of the same sign. In this case, two global extrema of L(s) + L(—s) exist at so = ±s* 
with s* > 0, such that two other families of kink modes bifurcate from so = ±s*. 

When V = V\(x), (fro = tanhx, and go = 1> the function L(s) is computed in the implicit form 



Clearly L(0) < and, as can be seen from Fig. ^(top left panel), L"(0) > for any Additionally, 
Fig. ^ (middle and bottom left panels) suggests that M(s) < and M(s) — > as s — > oo for any k. 
Therefore, there is only one kink mode that bifurcates from sq = 0, where V\{x) has a minimum. 

When V = V2(x), the function L(s) is computed in the implicit form: 



where k E M+. The above L(s) can be expressed as a generalized hyper geometric function, however, 
we will not reproduce the resulting expression here. Instead, we note that L(0) > and 



When k = 0, L"(0) = 2. By using the Laplace method for computations of the integrals, one can find 
that L"(0) = —4k -3 + 0(k -5 ) as k — > oo. Therefore, there exists Kq G R+ with Kq < Kg, such that 
L"(0) > for < k < Kq and L"(0) < for n > Kq . As can be seen from Fig. ^ (top right panel), 
K o = K o ~ 3.21. Additionally, Fig. ^ (middle and bottom right panels) suggest that no other extremal 
points of M{s) exist on s £ R for any k G R+. When < k < Kq, there exist three kink modes: two 
modes with sq = ±s* are associated with the global maxima of the effective potential M(s), while the 
mode with so = is associated with the local minimum of M{s). When k > Ko, no kink modes with 
so = is* exist but the mode at so = corresponds to the global maximum of the effective potential 
M(s). Hence, the structure of kink modes corresponds to a subcritical pitchfork bifurcation in the 
parameter k, such that three solutions exist (at so = and so = ±s*) for < k < and only one 
solution persists for k > kq. We point out that the effective potential M(s) gives a different prediction 
in comparison with the true potential V^ic) which possesses a minimum at x = and two maxima at 
x = ±2/k, for all k. 
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Figure 1: The effective potential evaluated numerically for V\{x) (left) and V^x) (right): the quantity 
L"(0) versus k (top panels) and the function M(s) for k = 1 (middle panels) and k = 4 (bottom 
panels). 
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3 Stability analysis of black solitons 



We first consider the spectral stability of black solitons in the NLS equation We will obtain 

conditions for spectral stability and instability of kinks and bubbles and then extend these conditions 
to kink modes of the GP equation (|1.2|) , In the end of this section, we will apply these conditions to 
the kink modes related to the two potentials Q1.4|l . 



Definition 3.1 Let 4>o( x ) be the black soliton of Definition \2. 5l The black soliton is said to be spectrally 
unstable in the time evolution of the NLS equation M.l\) if there exists an eigenvector (u, w) G L 2 (K, C 2 ) 
of the spectral problem 

L + u = —Xw, L-w = Xu, (3-1) 
for the eigenvalue X with Re(A) > 0, where 

L + = -\d 2 x + f{<&) - /(g ) + 2$/'($), L_ = -\d 2 x + /(0 2 ) - /(g ). (3.2) 
Otherwise, the black soliton is said to be spectrally stable. 



Remark 3.2 The spectral problem Q3.1|) arises in the linearization of the NLS equation (jl.lj) by using 
the ansatz 

u(x, t) = e ~ if( - qo)t <j> (x) + e xt [u(x) + iw(x)} + e xt [u(x) + iw(x)] 

It will be clear from analysis of the system Q3.1JI that the spectral instability of black solitons is always 
associated with a real positive eigenvalue A, while the spectral stability of black solitons (under a non- 
degeneracy constraint) corresponds to the case when the black soliton is a ground state of an equivalent 
variational principle. It is relatively straightforward to develop the nonlinear analysis for these two 
cases and to show that the spectral instability and stability of black solitons (under a non-degeneracy 
constraint) correspond to their orbital instability and stability. See [U I28| I31j for developments of this 
nonlinear analysis. 



Proposition 3.3 Let u(x) be a suitable function such that u(x) = ^/qoe tei± exponentially fast as 
x — > dboo. Let the renormalized energy E r [u] and momentum P r [u] of the NLS equation M-l\) be 
defined by 



E r [u] 



\u~r + 2 



</0 



M 2 



(f(qo)-f(q))dq 



dx, 



p rbA = 7; (uu x - uu x ) [ 1 - -j-^ ) dx. 
2/r V M J 



(3.3) 
(3.4) 



The family of dark solitons of Definition \2.1\ is a critical point of the Lyapunov functional A[u] 
E r [u] +vP r [u\. 



Proof. By direct differentiation, if u = U(x) satisfies the second-order ODE (|2.3|) with uj = f(qo), 
then the variational derivative £'^[n]| M= [/ + vP' r [u\\ u= u is zero. The second term in P r [u] is a Casimir 
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functional, which has zero variational derivative. Therefore, the Lyapunov functional A [it] is defined 
with accuracy to an arbitrary constant C in 



A[u] = E r [u] + vP r [u] + CS[u] , (3.5) 

where 

S[u] =~f C^-^)dx = [axg(u)] |™oo (3-6) 

is the total phase shift of u{x) on x G M subject to the non-zero boundary conditions on u{x) at infinity. 
In order to define the constant C uniquely, we add a constraint on the variational problem by requiring 
that if u = (fie tS , 6' = v(l — qo/cj) 2 ), and <f> = &(x) satisfies the second-order ODE (|2.7|) . then the first 
variation A'[0]|^ = $ is zero, where A[0] = A[it]| u=( x e ie, 0' = v(l — qo/(f> 2 ). By direct differentiation, it 
follows immediately that C = under this requirement. □ 



Remark 3.4 The renormalized momentum 1)3. 4|) was constructed in [22] as a difference between the 
standard momentum P[u] associated with a solution u(x) and the value P[«o] evaluated at the back- 
ground solution no = y / go"e* slgn ^^' So/ ' 2 , where the value Sq = S[u] is related to the total phase shift of the 
solution u(x). It was shown in ^[2B] that the renormalized momentum P r [u] computed at the family 
of dark solitons U(x) of Definition 12.11 defines the spectral stability and instability of dark solitons 
in the sense that the dark soliton is spectrally stable if P{.(v) > and unstable if P' r {v) < where 
P r (v) = P r [U]. (The degenerate case P' r (v) = corresponds to the dark solitons which are spectrally 
stable and orbitally unstable. Under the non-degeneracy constraint P' r {v) ^ 0, the spectral stability 
and instability corresponds to the orbital stability and instability, see |28j . In what follows, we will 
consider dark solitons under the non-degeneracy constraint P^(v) ^ Vi> S (— c, c).) We will show that 
the limit v [ is well defined and the quantity P/luJ.0 determines spectral stability and instability of 
kinks of Definition 12.51 The latter point is missed in the recent paper (SJ; where stability of kinks is 
considered. 



Lemma 3.5 Let U(x) be the family of dark solitons of Definition \2.1\ and f(q) be C 2 (M, + ). Then, 
(i) the function P r (v) = P r [U] is C 1 on v G (— c, 0) U (0,c) and (ii) the limiting quantity P^kio is 
well-defined. 



Proof, (i) By construction of dark solitons in Theorem 12.31 the function P r (y) = P r [U] is represented 
by 



i 



P r (v) = i/ (UU'-UU')(l-^)dx 



-v f $ 2 (x) ( 1 



2 

dx = vN(v) +q S(v), (3.7) 



<£ 2 (x 

where N(v) and S(v) is the total power and phase shift of the dark solitons: 

N(v)= [ (q Q -<l> 2 {x))dx, S(v)= [ 9'(x)dx = 0+-e-. (3.8) 
Jr Jr 

By the ODE theory for the system ((23|) - (|2~7|) with > on x G R, the map v >-* O) is C 1 on 

v G (— c, 0) U (0,c), such that N(v) and S(v) are smooth functions and so is P r (v). 
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(ii) We will show that the functions N(v), S(v) and P r {v) remains smooth in the limit v I 0. Let 
U(x) be a dark soliton for v G (0, c) according to Definition 12.11 and (j)o(x) be a black soliton according 
to Definition 12.51 Let us consider 

U(x) = U ^ - = U r (x) + iUAx), «€(0,c), 

v 

where U r {x) and ^(a;) are real- valued functions on x G R. By the construction of U{x) and </>o(a;) in 
Theorems 12.31 and 12. 6| it is clear that U r ,Ui G L°°(R) are continuous in v for v G (0, c). We need to 
prove that these functions remain continuous in v as w 1 0. By separating the real and imaginary parts 
in the ODEs (|2.3|) and ([2.10^ . we obtain the equivalent ODE system for U r (x) and Ui(x): 

vU{ + -0o (/($) ~ f{\U\ 2 )) + ^ + U r (/(go) - f(\U\ 2 )) = 0, (3.9) 
-0' o - vU' r + hj? + Ui (f(q ) - f(\U\ 2 )) = 0, (3.10) 

where |LT| 2 = (0 O + vU T ) 2 + v 2 U 2 . Let us rewrite the ODE (13111 as an inhomogeneous problem: 

(L- + L-)Ui = F-, (3.11) 
where operator L_ is defined by (|3.2jl and 

L- = f{\U\ 2 ) ~ /($), F- = -0o " vU'r- 

Since |C/| 2 (x), 4>q(x) approach go and U r (x), Ui(x), <po(x) approach some constants exponentially fast 
as \x\ — ► oo, it is clear that L_ is a relatively compact perturbation to L_ and L—,F~ G L 2 (R) for 
v G (0, c). By continuity of the solution <po(x) = lim £7 (a;) in Theorem l2.61 we know that ||L_ ||x,oo = o(l) 

and u||C/^||lcx) = o(l) as v j 0. 

Since L_0 O = and O G L°°(R), then V/ G L 2 (R) there exists Li 1 / G L°°(R) if and only if 
(00,/) = 0- The following computation shows that this condition for / = F- — L_[/j is equivalent to 
the ODE ()3.1U|) and is thus satisfied on v G (0, c): 

-(0o, 0' o ) - v(4> , U' r ) + (0o(/(0o) " /(<Zo), Ui) + (0o, ^ (/(g ) - /(|t/| 2 ))) 
= -(0o, 0' o ) " «(fo,Eft + i^,^) + (0o, ^ (/(go) - /(l^| 2 ))) 
= (0o,-0o-< + if7r + C/ i (/(go)-/(|C/| 2 ))) =0. 

Therefore, the operator LZ 1 '■ L 2 (R) ^ L°°(R) is well-defined for the inhomogeneous problem ()3.11|) 
on v G (0, c), which can be rewritten as follows: 

Ui = - (L- (<Po + vufj . 

Since = o(l) and = o(l) as v j 0, there exists a solution Ui G L°°(M) uniformly in 

v G [0, c), such that C/i|«|o = —LZ 1 ^ G L°°(R). Therefore, the function lmU(x) is smooth as v 1 

and Im^C7(a;)|^o = Ui|„| = --ECVo- 
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We can now use the fact that i>||ft||L°° = 0(v) as v J. 0. Since /(g) is C 2 (M + ) and v ||ft||L°° = °( v ) 
as v I 0, there exists a function g(cpQ,vU r ) for small ut/ r such that 

/((0o + ^ft.) 2 ) - f{4>l) - 2vfoU r f(<j>l) = vU r g(cl> ,vU r ), 

where \\g{4>o,vU r )\\L°° = o(v) as v J. 0. Using these facts, we rewrite the ODE ()3.9|) as an inhomogeneous 
problem: 

(L+ + L+)ft = F+, (3.12) 
where operator L + is defined by (|3.2j) and 

L+ = /(|£/| 2 ) - f(<j>l) + fogifavUr), F + = vUi + (/((fo + vU r ) 2 ) - f(\U\ 2 )) . 

It is clear that L + ,F+ G L°°(R) for t> G (0,c). Since L+0' o = and $,(x) G L 2 (IR), then V/ G L°°(R) 
there exists f G L°°(R) if and only if ((f>' , /) = (by the Fredholm's Alternative). The following 
computation shows that this condition for / = F + — L + U r is equivalent to the ODE ()3.9j) and is thus 
satisfied on v G (0, c): 

v(<P' , %) + I (0^o, (/fog) - /(Ift 2 ))) + (V (/($) + 20g/'(0g) - /(|f/| 2 )) , ft) 
= t,(^, ft) + i (tfrfo, (/fog) - f(\U\ 2 ))) + ft-) + (0o, ft (/(ft) - /(Ift 2 ))) 

= (Vo, < + ~0 O (/($) - /(Ift 2 )) + ~U? + ft (f(q ) - f(\U\ 2 ))) = 0. 

Therefore, the operator L^ 1 : L°°(R) h- > L°°(R) is well-defined for the inhomogeneous problem ()3.12j) 
on v G (0, c), which can be rewritten as follows: 

ft = - (£+ + L+y 1 [vUi + i^o (/((& + vU T f) - /(|lf ))) ■ 

Since [|L + [|l<xj = o(v) and tt[|ft[|.z,°o = 0(v) as u 1 0, there exists a solution £/ r G L°°(R) uniformly 
in v G [0, c), such that ft|«j.o = (the homogeneous solution (p' (x) is removed from U r (x) due to the 
symmetry in U(x)). Therefore, the function KeU(x) is smooth as v [ and Ked v U(x)\ v io = ft|«j.o = 0. 

As a result, the map v \— ► U is C 1 on i> G [0, c), such that iV(v) and S(v) are smooth functions as 
v I and so is P r (v). □ 

Corollary 3.6 XYie following identities are true for v G (— c, 0) U (0, c) 

= i [ (U'd v U - U'd v U) dx = 2 [ (Ret/' lmd v U - JmU' Red v U) dx (3.13) 

and 

PXlo = 2 {<P' ,Imd v U\ vl0 ) = N\ vi0 + q S'\ vio . (3.14) 
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Proof. By Lemma 13.51 the quantity P' r {v) is continuous on v £ [0, c). The first identity in (|3.13|) follows 
by direct differentiation: 



p' r (v) = 1 f (u'd v U + Ud v U' -U'd v U - Ud v U') dx + ^ [ d v 



ago 
2 ' 



17' 








f7 






§) 


■( 





X— > — oo 



(U'd v U -U'd v U) dx. 



Other identities follow by the substitution U{x) = Rei7(x) + iImU(x) } by the smoothness of U(x) with 
respect to v £ [0, c) and by the relation (|3.7|) . □ 



Example 3.7 Following Example 12. 15( we consider the cubic NLS with f(s) = s and qo = 1. By using 
the exact solution 1)2.17(1 . we find for v £ [0, 1) 



N(v) = 2VT - v 2 , S(v) = -2arct 

such that 

p' r ( v ) =4Vl-v 2 , S'(v) 



an- 



and PXio = 4, S"|^o = 2. 



Lemma 3.8 (%) Let (j)o(x) be the kink of Definition \2. 3i Then, the spectrum of L + in L 2 (M) consists 
of the positive continuous spectrum bounded away from zero by 2c 2 , the kernel with the eigen] 'unction 
<p'o(x) and, possibly, a finite number of positive eigenvalues in (0,2c 2 ). The spectrum of L- in L 2 (M) 
consists of the non-negative continuous spectrum and a single negative eigenvalue. 

(ii) Let 4>o(x) be the bubble of Definition 1^751 Then, the spectrum of L + in L 2 (M) consists of the 
positive continuous spectrum bounded away from zero by 2c 2 , the kernel with the eigenfunction 4>q(x), 
a single negative eigenvalue and, possibly, a finite number of positive eigenvalues in (0,2c 2 ). The 
spectrum of L_ in L 2 (]R) consists of the non-negative continuous spectrum. 



Proof. Since 4>q{x) — > qo exponentially fast as |x| — > oo, L± in (|3.2|) are self-adjoint Schrodinger 
operators on the domain i7 2 (R) C L 2 (IR), which have absolutely continuous spectrum a c (L±), a finite 
number of isolated eigenvalues of finite multiplicities o~ p (L±), and no embedded eigenvalues or residual 
spectrum [Tl]. By the Weyl's Essential Spectrum Lemma, a c (L + ) > 2go/'(?o) = 2c 2 > and <t c (L_) > 
0, such that the continuous spectrum of L + is bounded away from zero and the continuous spectrum 
of L_ touches zero. Moreover, L + cjJ (x) = and L_0o(x) = due to the translational and gauge 
symmetries of the NLS equation (|1.1|) . such that L+ has a simple kernel in L 2 (R) while L_ has no 
kernel in L 2 (M). 

(i) In the case of kinks, (j)o( x ) has a single zero on x £ M. By the Sturm Nodal Theorem, a p (L + ) 
contains no negative eigenvalues and o" p (7_) contains exactly one negative eigenvalue. 

(ii) In the case of bubbles, 4>o{x) has no zeros on x £ R. By the Sturm Nodal Theorem, a p {L + ) 
contains exactly one negative eigenvalue and cr p (7_) contains no negative eigenvalues. □ 
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Lemma 3.9 Consider the constrained space 

X c = {w£ H 2 (R) : (w, <f/ ) = 0} , (3.15) 

where (f>o{x) is a black soliton of Definition \2.b\ In the case of kinks, the operator L_ has exactly one 
negative eigenvalue in X c i/P/l^o < and no negative eigenvalues if P^| v j.o > 0, where Pr\vio is defined 
by Lemma \'J. 51 In the case of bubbles, the operator L_ is non-negative in X c . 

Proof. We consider a constrained variational problem: 

(£_ - fj,)w = -u(f)' , w E H 2 (R), fx <£ <t(L_), (3.16) 

where [i is the spectral parameter and v is the Lagrange multiplier. By Lemma 13.81 there exists a 
unique solution w E H 2 {M) for any \x E (/^OjO), where [xq is the only negative eigenvalue of L_ in the 
case of kinks and /xo = — oo in the case of bubbles. By the standard variational theory (e.g. see [T|l3T]). 
the smooth function g(/x) = — (4>' , (L_ — //) _1 <^>q) is decreasing on fj, E (fJ.Q,0) from lim g(fi) = +oo in 

the case of kinks or lim g(p) = in the case of bubbles. Therefore, operator L_ is non-negative in 

X c in the case of bubbles. In the case of kinks, operator L_ is non-negative in X c if lim g(p) > and 

it has a single zero on fi E (fJ-o,0) if hmg(/j,) < 0. We need to show that limg(fi) = P' r \ v iQ. 

/40 //To 

It is clear from the solution of the variational problem (|3.16|) that w E H 2 (R) C L°°(R) uniformly 
in [A E (a*0)0]. Since (<^o,$)) = 0> there exists a solution of the inhomogeneous problem L_wq = — (f>' 
in wq E L°°(M) and, moreover, it follows from smoothness by Lemma 13.51 that 

= Ui\ vlQ = lmd v U(x)\ vi0 + c(j} , 

where c is arbitrary. By the relation (|3.14|) . we obtain lim g(f/,) = (<j)' ,lmd v U\ v io) = nP^Uio- ^ 

Remark 3.10 It is interesting to note that the inhomogeneous equation L_wo = — <Pq admits another 
solution wq = x4>o(x) ^ L°°(R), which may lead to miscalculation 

limg(n) + {(f>' ,x<p ) = -N r \ v ± > 0. 
/tjO z 

The linearly growing solution x(f>o(x) ^ L°°(R) is however irrelevant for the variational problem (J3.16|) . 

Theorem 3.11 (i) Let 4>q{x) be the kink of Definition \2.3[ Then, it is spectrally stable if P{.\vio > 
and unstable if -P^lv-io < 0. In the latter case, the spectral problem hS. 1\) has exactly one real positive 
eigenvalue A. 

(ii) Let <f>o(x) be the bubble of Definition \2.5\ Then, it is spectrally unstable with exactly one real 
positive eigenvalue A in the spectral problem 

Proof. Let A be a non-zero eigenvalue of the spectral problem (|3.1j) with (u, w) E H 2 (R, C 2 ). Then, L+ 
is invertible in X c defined by (|3.15|) and the component w(x) E X c can be found from the generalized 
eigenvalue problem 

L_w = -yLl 1 w, 7 = -A 2 , weX c . (3.17) 
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Due to the equivalence above, all non-zero eigenvalues A of the spectral problem (|3,1|) can be recovered 
from the non-zero eigenvalues 7 of the generalized eigenvalue problem (J3.17|) • The operators L± satisfy 
properties P1-P2 of the recent paper JSJ. Even though L_ has no spectral gap near the origin, one can 
shift the generalized eigenvalue problem to the equivalent form, 

(L_ + SL+^w = (7 + <5)L+V < S < <5 , (3.18) 

where 5q is the distance from 7 = to the first negative eigenvalue 7 if it exists or 5q = 00 if not. Now 

a c (L_ + 5L- 1 ) > ^2 > 0, 

and the operator L_ = L_ + has the spectral gap near the origin. 

(i) In the case of kinks, the operator L^ 1 is positive in X c . By Theorem 3 of jH], the problem ()3.18j) 
has no eigenvalues 7 G C with Im(7) 7^ 0, has no eigenvalues 7 G R such that (w, L7}w) < and 
has exactly ./V = dim(ff7 eigenvalues 7 < 0, where H~ 1 C X c is the invariant negative 

Li— -f-Oi/_|_ i/— -\~OLi\ 

subspace of X c with respect to L_ + ^L^ 1 . We will prove that N = dim(H L +gL -i) = dhn(H£_) 

for sufficiently small 5 > 0. Indeed, continuity of isolated negative eigenvalues of L_ in 5 follows 
by the perturbation theory since 5L^_ is a relatively compact perturbation to L_ in X c . Therefore, 
dim(H~ +SL -i) > dim(iJ^ ) for sufficiently small 5 > 0. Consider a splitting X c = ^7^_ © H^_, 

where H£ (Ht ) is negative (non-positive) invariant subspace of X c with respect to L_, such that 
dim(H^ ) < 00 and dim(H~£ ) = 00. Since L^ 1 is strictly positive operator in X c , we have 

V<5 > 0, Vu> G F+_ : (u>, (L_ + taT. 1 )™) > £(«>, L; 1 ^) > 0. 

Therefore, the operator (L_ + 5LJ1 1 ) is strictly positive on H£ and dim(i?~ < dim(ii^ ). We 

have thus proved that dim(//^ = dim(.ff~£ ). By Lemma 1331 dim(H^ ) = 1 if P^\ v iq < and 

dim(Hj^ ) = if > 0. In the former case, N = 1 and the kink is spectrally unstable. In the 

latter case, N = and the kink is spectrally stable. 

(ii) In the case of bubbles, the operator L_ has no negative eigenvalues in X c by Lemma 13.91 
Therefore, Theorem 3 of (§] guarantees that the generalized eigenvalue problem ()3.18j) has no eigenvalues 
7 G C with Im(7) 7^ 0, has no eigenvalues 7 G M+ such that (w,L~Z w) < and has exactly N = 
dim(H~_ 1 ) eigenvalues 7 G R_ with (wjL^w) < 0, where H~_ 1 C X c is the invariant negative 

subspace of L'Z 1 . Since all eigenvectors of L + for non-zero eigenvalues are orthogonal to (p' Q and belong 

to X c , it follows immediately that N = dim(H~_ 1 ) = 1. □ 

L + 

Remark 3.12 The statement (i) of Theorem 13.111 extends the stability-instability theorem in |2H] 
from dark solitons with v 7^ to kinks with v = 0. This connection is missed in the recent paper 
where the same result was obtained by the Vakhitov-Kolokolov method (used in pQ) and the variational 
principle (used in [I]) without the proof that limg(^) = ^Pr\vio- 

The statement (ii) of Theorem 13.111 was proved differently in 0] by using a variational technique. We 
note that P r {v) < for v G (0, c) and limP r (v) = in the case of bubbles of Definition 12.51 Therefore, 

if P/|uio 0) then -P^Uio < 0) such that the statement (i) for kinks extends formally to the statement 
(ii) for bubbles. 
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Definition 3.13 Let <$> e {x) be the kink mode of Definition ^. SI The kink mode is said to be spectrally 
unstable in the time evolution of the GP equation hi. 2^ if there exists an eigenvector (u,w) E L 2 (W, C 2 ) 
of the spectral problem 

C + u = —Xw, C-W = An, (3.19) 
for the eigenvalue X with Re(A) > 0, where 

= ~\dl + M) ~ /(go) + Wlfm + eV(x), £- = ~\dl + f(<%) - /(go) + eV{x). (3.20) 
Otherwise, the kink mode is said to be spectrally stable. 

Theorem 3.14 Let 4> e (x) be the kink mode of Definition \2.tft Assume that the operators C± have 
n± negative eigenvalues and empty kernels in L 2 (IR). Assume that all embedded (purely imaginary) 
eigenvalues of the spectral problem hS.iy^ are algebraically simple. Then, the spectral problem \3.1fJ\) 
for (u, w) G L 2 (M, C 2 ) has exactly N c complex eigenvalues A in the first quadrant, N~ purely imaginary 
eigenvalues A with Im(A) > and (w,£7 w) < 0, and N r = iV+ + N~ real positive eigenvalues A, 
where N^~ corresponds to eigenvalues with (tu,£7 w) < and N^~ corresponds to eigenvalues with 
(w,£^_ w) > 0, such that 

N+ + Nf +N c = n+, N r 7 + Nr + N c = n_ , (3.21) 

where multiple eigenvalues are accounted up to their algebraic multiplicities. 



Proof. Since V(x) — > and <t> 2 {x) — > qo exponentially fast as \x\ — ► oo, operators C± have the absolutely 
continuous spectrum such that a c {C + ) > 2c 2 > and cr c (£_) > 0. Since the kernel of £+ is empty in 
L 2 (IR) by assumption, the generalized eigenvalue problem ()3.17j) for operators C± is rewritten in the 
unconstrained space: 

£-w = "/£r 1 w, 7 = -A 2 , w£H 2 (R). (3.22) 

The kernel of £_ is empty in L 2 (M) since C-(p t = and e G L°°(K), <j) t <£ L 2 (R). The generalized 
eigenvalue problem (|3.22j) can be rewritten in the equivalent form, 

{C- + 5Cr 1 )w = {~/ + 5)Cr 1 w , (3.23) 

where 5 > is sufficiently small. Properties P1-P2 of jHj are satisfied and Theorem 3 of gives the 
relations 

N+ + Nr + N c = dim(#- ! ) , iV r 7 + Nr + N c = dim(i? 7 ) 

for sufficiently small 5 > 0, where H~i and iJ7 a are invariant negative subspaces of L 2 (R) with 

respect to £7 and £_ +&C77 1 respectively. It follows immediately that dim(H~ 1 ) = n + . By continuity 

of eigenvalues and the relative compactness of Cr 1 with respect to it follows that dim(ff7 ) < 

dim(i/7 We shall prove that dim(i/7 xr -i) = dim(i?7_). The operator £_ + 5£r 1 may have 

additional negative eigenvalues compared to operator £_ if and only if some eigenvalues bifurcate as 
5^0 from the end point of the continuous spectrum of £_ by means of the edge bifurcation [1111201 US] . 
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In order to analyze the edge bifurcation, we rewrite the eigenvalue problem (£_ + 5C + )w = (jlw in the 
equivalent form: 

(£ + SM) w = iiw, wa 2 (l), 

where 

C = C^+5(2c 2 -\dl S j \ M=£+^f(tf)-f(q ) + 2<p 2 J>(tf) + eV(x))(2c 2 -^ * 

where M. is a relatively compact perturbation to the unbounded operator C. The continuous spectrum 
of C is bounded from below by cr c (£) > ^L-. By the theory of edge bifurcations (see review in [10]), 
the new eigenvalue \x = \xs, if it bifurcates from the end point of cr c (C), has the expansion fj,$ = 
ipx — a5 2 + 0(<5 3 ), where a is positive constant. Therefore, there exists sufficiently small 6 > 0, such 
that fj>s > 0. As a result, the edge bifurcation does not change the number of negative eigenvalues of 
£_ + SC^ 1 compared to £_ and dim(H~ +sc -i) = dim(iT^ ) = n_. □ 



Theorem 3.15 Let (j>o(x) be the kink of Definition s.^ and M'(s) be defined by \2.13\) . Let 4> e (x) be 
a unique continuation of the kink 4>o(x — sq) in Theorem \2.12\ from the root sq such that M'(sq) = 
and M"{sq) ^ 0. Then, the operators C± have n± negative eigenvalues and empty kernels in L 2 (W) 
for sufficiently small e with n+ = 1, n_ = 1 for M"(sq) > and n + = 0, n_ = 1 for M"{sq) < 0. 



Proof. It follows by Theorem 12.121 that (j) e = 4>q(x — s) + e<pi(x) + <p(x, e, s) and s = sq + s(e), where 
= o(e) and \s\ = o(l) as e — > 0. Therefore, operators C± in (|3.2()j) are represented by 

C± = L± + eM± + M±, 

where L± are given by ()3.2|) . M± are given by 

M+ = V(x) + 6foVif(4) + 4$<Pif" fao), M- = V{x) + 20oVi/(<Ao)> (3-24) 

and [|M±||lcc = o(e) as e 0. We note that eM + + M+ = eV(a?) + D^N^s^e), where ZA^iV is 
the Jacobian of the nonlinear function in (j2.14j) . By using the inhomogeneous equation (j2.16j) for the 
correction term cpi(x), we compute 

(</>' (x - s ), M + ^' Q (x - s )) = -(#,, Vfa) ~ (<4 L+iffj) = -\m"{s ). (3.25) 

By the regular perturbation theory, the zero eigenvalue of L + becomes a non-zero eigenvalue X t of £+ 
for small e, such that 

M"(s ) 



A e = eAi + A, Ai 



where A = o(e) as e — > 0. Since the zero eigenvalue of L+ is simple and positive eigenvalues of L + are 
bounded away from zero, the kernel of is empty, such that n+ = 1 for M"(sq) > and n+ = for 
M"(so) < for sufficiently small e. By the Implicit Function Theorem applied to (j>e(x) (since (f>o(x) 
has only one simple zero at x = 0, the function 4> £ (x) has only one node oniGl for sufficiently small 
e. We recall that = and <p t £ L°°(R), <p t <£ L 2 (R). By the Sturm Nodal Theorem, the kernel of 

£_ is empty and n_ = 1 for sufficiently small e. □ 
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Remark 3.16 Let / e C 2 (R + ). Then \\M±\\ L oo = 0(e 2 ) and A = 0(e 2 ) as e -> 0. 

Corollary 3.17 The kink mode with M"(sq) < is spectrally unstable with exactly one real positive 
eigenvalue A in the spectral problem \3.iy\) for sufficiently small e. The kink mode with M"{sq) > 
may have up to two unstable eigenvalues A in the spectral problem h'J.iy\) . 

Proof. If M"(sq) < 0, then n + = 0, n_ = 1 and the count of eigenvalues ()3.21j) gives iV+ = = 
N c = and N~ = 1. If M"{sq) > 0, then n + = ra_ = 1 and the count of eigenvalues may give 
either Nf + N c = 1, N+ = N~ = or Nr = N c = 0, = A^ T 7 = 1. In the cases N c = 1 or 
= iV r 7 = 1, there are two unstable and no embedded eigenvalues in the spectral problem (|3.19j) . In 
the case iV~ = 1, the pair of embedded eigenvalues is simple, such that the last assumption of Theorem 
EH is satisfied. □ 

Remark 3.18 Asymptotic approximations of eigenvalues A and precise statements on unstable eigen- 
values in the case M"(so) > are obtained in Section 4 under non-degeneracy assumptions i^k-io 
and S'lvio ^ 0. By Corollary 14.121 the case Nf = N c = 0, = N~ = 1 occurs for P' r \ v is < and 
the case N+ = N~ = Nr = 0, N c = 1 occurs for P^lwio > 0. 

Example 3.19 Continuing Examples 12 . 1 51 and 13 . 71 we consider the cubic NLS equation with /(s) = s, 
qo = 1 and P' r \ v io = 4 > 0. When the potential V{x) is even with V(—x) = V(x), the function M'(s) 
is represented by M'(s) = L'(s) — L'(—s), where L(s) is given by (|2,18|h One family of kink modes 
with s = always bifurcates with M"(0) = 2L"(0). 

When V = Vi(x), it follows from Fig. □ (top left panel) that L"(0) > for any k / 0. The kink 
mode with sq = corresponds to the minimum of M(s) and it is unstable with two complex conjugate 
eigenvalues, according to Remark 13.181 

When V = V^x), it follows from Fig. ^ (top right panel) that there exists < /-to < °° such that 
L"(0) > for < k < kq and L"(0) < for n > kq. For < k < Kq, a pair of kink modes bifurcates 
from so = ±s* with M"(sq) < 0. These modes correspond to the maxima of the effective potential 
M(s) and they are unstable with one real eigenvalue, according to Corollary 13.171 The kink mode 
with so = corresponds to the minimum of the effective potential M(s) and it is unstable with two 
complex eigenvalues for < k < Kq. On the other hand, it corresponds to the maximum of M(s) and 
it is unstable with a simple real positive eigenvalue for k > kq. This scenario indicates the subcritical 
pitchfork bifurcation at k = kq . We will illustrate this bifurcation in Section El 

4 Eigenfunctions and eigenvalues of kinks 

We report here asymptotic analysis of the spectral problem ([3.19ft in the limit of small A and e. This 
asymptotic analysis is needed to complete the stability analysis of kink modes with M"(sq) > which 
is not conclusive in Corollary 13.171 We will show that if the kink is stable in the linear problem 1)3 
for e = 0, then the pair of zero eigenvalues of the spectral problem ()3.19j) at e = splits into a pair 
of purely imaginary eigenvalues at O(e) and bifurcates into a quartet of four complex eigenvalues (two 
of which are unstable) at 0(e 3 / 2 ). These eigenvalues A for small e correspond to the eigenvectors 
(u,w) E L 2 (R, C 2 ), persistence of which in e follows by Theorem 13.141 
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From a technical point of view, our analysis is complicated by the fact that the eigenvalues A £ iM. are 
embedded into the continuous spectrum of the non-self-adjoint problem H3.19JI . Since we are interested 
in specific information about bifurcations of eigenvalues near the point A = 0, we have to abandon the 
generalized eigenvalue problem (|3.22|) and to work directly with the non-self- adjoint eigenvalue problem 



One way to deal with this problem is to introduce exponential weights which move branches of the 
continuous spectrum from the imaginary axis (see [3E1 and references therein). However, there are two 
branches of the continuous spectrum, and, independently of the weight parameter, one branch moves to 
the left and the other branch moves to the right of the imaginary axis. Any eigenvalues that bifurcate 
off the imaginary axis may become resonant poles when the weight parameter is sent to zero unless 
specific information about the decay rate of eigenfunctions is available. However, if this information 
were available, one could avoid the technique of exponential weights and perform a direct analysis of 
the eigenfunctions and eigenvalues of the problem (|3.19|) . 

Another way to deal with this problem is to consider the Evans function with careful analysis of 
fast and slow decaying solutions (see [HI and reference therein). By using the Gap Lemma, the Evans 
function can be appropriately extended across the continuous spectrum with a full account of the branch 
points on the imaginary axis. Information about small eigenvalues is drawn from the derivatives of 
the Evans function with respect to A and e near A = and e = 0. However, computational formulas 
become more and more involved when higher-order derivatives of the Evans function are needed and 
this has caused some miscalculations in the past (see Appendix A in |44j ) . 

Our treatment of the problem brings together the analysis of fast and slow decaying solutions in 
the two approaches above. To avoid complications, it is based on direct analysis of eigenfunctions and 
eigenvalues of the spectral problem (|3.19j) expanded in powers of e 1 / 2 . We will obtain a characteristic 
equation for small eigenvalue A versus small parameter e. 

Lemma 4.1 Let 4>q(x) be the kink of Definition \2.b\ and U{x) be the dark soliton of Definition \2.i\ for 
v > 0. Let operators L± be defined by US. 2\) . The uncoupled homogeneous problems 



admit four linearly independent solutions: 

(i) exponentially decaying eigenfunction uq = (Pq(x) in L 2 (M) 

(ii) bounded eigenfunction wq = 4>q(x) in L°°(M) 

(Hi) unbounded linearly growing solution wo = x<fio(x) — Imd v U(x)\ v iQ 

(iv) and an unbounded exponentially growing solution uq 

The uncoupled inhomogeneous problems 



L + u = 0, 



L-Wq = 



L+ui = -w , 



L_w\ = uq 



admit solutions in the same order: 



(i) bounded eigenfunction w± 



lmd v U(x)\ vl0 in L°°(R) 
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(ii) a bounded eigenfunction U\ in L c 

(iii) an unbounded exponentially growing solution u\ if S'\ v \o 7^ 

(iv) and an unbounded exponentially growing solution w± 

The uncoupled inhomogeneous problems 

L + u 2 = —wi, L-W2 = u\ 
admit no solutions in L°°(R) if P^\ v \0 7^ 0. 

Proof. It follows from the proofs of Lemmas 13.81 and 13.91 that 

L + (p' = 0, L_0 O = O, L^xcpo = -(p' , L-Imd v U(x)\ v i = -<f>' , (4.1) 

which proves (i)-(iii) for uq and wq and (i) for w\. Existence of an exponentially growing solution 
uq in (iv) follows from the fact that the Wronskian determinant of two linearly independent solutions 
of L + uq = is constant in x. Existence of bounded solution u\ in (ii) follows from the fact that 
(4>' Q , 0o ) = 0, such that LZ^4>q 6 L°°(R). The solution u\ in (iii) grows exponentially since the Fredholm 
Alternative is not satisfied for wq in (iii): 

{^xMx) ~ Imd v U(x)\ vl0 ) = ±N\ vl0 - -PXio = -jS% l0 + 0, (4.2) 

where the relation (j3.14j) has been used. Existence of an exponentially growing solution w\ in (iv) 
follows from the fact that LZ 1 ^ has the same exponential growth in x as uq in (iv). The solution u 2 
of L + U2 = ImdyU ~{x)\vio grows exponentially since the Fredholm Alternative is not satisfied: 

(4>' ,lmd v U(x)\ v io) = -PXio + 
The solution w 2 of L_w 2 = u\ in (ii) grows linearly due to the same reason since 

(00, ui) = ~(L + ui,ui) 7^ 0. 

The last inequality is due to the non-negativity of L + for kinks and the orthogonality of the odd 
function u\ to the even function 4>' of the kernel of L + . □ 

Definition 4.2 Let A be fixed in the strip {A G C : < ReA < c 2 } and define k±{\) from the roots of 
the characteristic equations 



ReK± > : k 2 ± = 2c 2 ( 1 ± \l 1 - ^ ) , (4.3) 



such that k + k„ = 2A and k + = ^ /4c 2 



Remark 4.3 The roots k± can be expanded in the Taylor series near A = 0, such that 



k+(X) = 2c ( 1 - + 0(A 4 ) , «_(A) = - 1 + -j + 0(A 4 ) ) . (4.4) 
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Lemma 4.4 Let 4> e (x) be the kink mode of Definition ^. (A in the domain —eg < e < eo for some cq > 0. 
There exist four fundamental solutions (u, w) of the spectral problem \3.1fA) for any A in the strip 
{A G C : < ReA < c 2 }, such that 



u \ I K± 



W J , \ —K T 

and 

' U \ _/ K± 
W J \ -K T 



e K±x as x -> -co (4.5) 
e~ K±x as x -> +oo (4.6) 



Proof. For any A in the strip {A G C : < ReA < c 2 }, the two roots K+(A) and «_(A) of the 
characteristic equations (|4.3j) are distinct and Ren± > 0. Existence of four linearly independent 
solutions with the exponential tails in Q4.5|) and (|4.6|) follows by the Coddington-Levinson's Theorem 
for ODEs ^U] under the condition that V{x) — ► and 4> 2 {x) — > qo exponentially fast as |x| — ► oo. □ 



Definition 4.5 T/ie determinant of the four fundamental solutions in Lemma \4-4\ Q t an V x 6 
called the Evans function E(X, e) o/ the spectral problem AS.lty . namely 



E(X, e) = det 



u + 


u + 


u_ 


U- 


u> + 


u' + 


u'_ 


u'_ 


w + 


w + 


W- 




w> 


w> 


w'_ 


w'_ 



(4.7) 



Remark 4.6 Because the Wronskian determinant of any four particular solutions of the ODE (|3.19j) 
is independent of x, the values of E(\, e) are independent of x. 



Lemma 4.7 Let <j) e {x) be the kink mode of Definition ^. fl| A = ^K-K + and k + = \j^c 2 — k 2 . The four 

fundamental solutions and the Evans function E(X, e) of the spectral problem \3.iy\) are analytically 
continued in variable k_ near k_ = for sufficiently small e. 

Proof. Let us unfold the branch point A = with the transformation 

A = c 2 = <±£. (4.8) 

The spectral problem (|3.19j) is rewritten explicitly as follows: 

[-d 2 x + k\ + k 2 _ +2V + {x)]u = -k + k_w, [-d 2 x + 2V-{x)} w = k + k_u, (4.9) 

where 

V+(x) = f(<f> 2 ) - f(q ) + 2<f> 2 f'(4> 2 ) ~ 2q f'(q ) + eV(x), V-(x) = f(4> 2 ) - f(q ) + eV(x). 

Since the ODE system (|4.9j) depends analytically on («+, k_) G C 2 and the boundary conditions (|4.5jl - 
(|4.6I) are also analytic in variables the four fundamental solutions are analytic on k_) G 
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C 2 and so is the Evans function E(X, e) as a determinant of analytic functions at any fixed The 
unfolding transformation (|4.8|) implies that the parameter c 6 C is arbitrary. Since c G M+ is fixed, this 

leads to the constraint n + = yj 4c 2 — k 2 , which is locally analytic near k_ = 0. The analytic functions 

on G C 2 with the locally analytic constraint k + = yj 4c 2 — k 2 are locally analytic functions 

in a neighborhood of k_ =0 for any sufficiently small e. □ 

Remark 4.8 The Evans function was constructed in ^Hj for the cubic NLS equation with a pertur- 
bation. It was discussed in ^5] under a general set of assumptions that the function E(X) is analytic 
in a small domain near A = with ReA > 0, its zeros coincide with eigenvalues A with the account of 
their algebraic multiplicities, and it is analytically continued in the variable k_ near the point A = 
= 0). Our analysis of Lemma 14.71 is different. It is based on the unfolding transformation (|4.8[) 
similarly to the recent work |12j . 



Example 4.9 Continuing Example l2.15l we compute the Evans function E(\) explicitly for the cubic 
NLS with f(s) = s and qo = 1. There exist explicit solutions of the spectral problem IjH.ljl for the cubic 
NLS (see, e.g. )■ By using these solutions, we obtain the explicit representation of the eigenvectors 
(u,w) satisfying the boundary conditions Ij4.5|) and (|4.6j) : 

u± = -gljT — e K±x ^sech 2 x + k ± tanhx - , w± = — e K±x (k± - 2 ta.nh x) (4.10) 

and 

u± = — -— — e~ K±x ( sech 2 x - k± tanhx - -K± J , w± = -77— ^ — e~ K±:E (k ± + 2tanhx) . 
2 + k± \ 2 J 2 + k± 

(4.11) 

The Evans function £/(A) is computed explicitly as the determinant of the four fundamental solutions 
in the form 

E(A) - (^ww - (412) 

such that -E(A) = 8A 3 (l — A + 0(A 2 )) as A — > with ReA > 0. The validity of all explicit formulas has 
been confirmed by means of the Wolfram's Mathematica. 



Theorem 4.10 Let f(q) beC°°(R + ) andV(x) beC 2 (R) satisfying 1731) . LetPXio + andM"(s ) = 
in Theorems \2.12\ and \3. 1 li Letj = T(e) be an eigenvalue of the generalized eigenvalue problem \3.22)) 
for sufficiently small e, such that w G L 2 (IR), Imr < 0, and limr(e) = 0. Then, the spectral problem 

\3.iy\) for sufficiently small e admits two eigenvalues A = ±A(e) = iy 7 — r(e) with (u,w) G L 2 (R, C 2 ), 
ReAImA > 0, and lim A(e) = 0, such that 

(i) A(e) is infinitely smooth with respect to e 1 / 2 . 

(ii) (u,w) is infinitely smooth with respect to e 1 / 2 and 

]xolu(x) = 4>b(x), limu;(a;) = 0, (4-13) 

up to an arbitrary multiplicative factor. 
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(Hi) (u,w) admits an asymptotic expansion as A — ► ; ReA > for large ±x 3> 1: 



a± (1 + 0(A)) ( X t + °S\ )(^~ + 0(Az)* ) , (4.14) 



c 



w ; ^ v ' VJJ V -2c + 0(A 2 ) 
where a± are some constants and A = A(e). 

Proof. Under the conditions of the theorem, there must exist a zero of the Evans function E(X, e) at 
A = A(e) for sufficiently small e, such that ReAImA > and lim A(e) = 0. 

(i) By Lemma 14,71 the Evans function E(X,e) is analytically continued in k_ = A/c + 0(A 3 ) near 
A = 0. It is also infinitely smooth in e near e = 0. (Indeed, the potential terms V±(x) in the 
representation (|4.9|) are infinitely smooth in e and exponentially decaying as \x\ — > oo.) In addition, 
the Evans function E(X, e) has the following properties: 

E(X, 0) = a\ 3 + 0(A 4 ), E(0,e)=0, 

where a is a numerical constant. According to Lemma l4,H the triple root of E(X, 0) corresponds to the 
generalized kernel (i) due to translational invariance and the bounded solutions (ii) due to the gauge 
invariance. The former subspace results in the double root A = of E(X,0) (as a proper eigenvalue), 
while the latter subspace results in a single root k_ = (A = 0) of E(X, 0) (as a proper resonance) |2*U] . 
By Lemma f4.ll a / is equivalent to the condition -P^kio ^ 0. The constraint E(Q, e) = follows 
from existence of <p e G L°°(R) such that C-(p e = 0. As a result, the Evans function is expanded near 
A = and e = as follows: 

E(X, e) = A (aA 2 + f3e + 0(A 3 , Ae, e 2 )) , (4.15) 

where (3 is another numerical constant. Since the gauge invariance is preserved while the translational 
invariance of dark solitons is destroyed by Theorem 13.151 (3 ^ is equivalent to the condition that 
M"(sq) 0. It follows from expansion (|4.15j) and the smoothness of E(X,e) that the root A(e) is 
infinitely smooth in e 1 / 2 . 

(ii) By the construction of the Evans function E(X, e), the eigenvector (u, w) is spanned by the four 
fundamental solutions in Lemma 14.41 By Lemma 14.71 these solutions are analytic in k_ near K- = 0. 
By Remark 14.31 k_ = A/c + 0(A 3 ) near A = and by (i) of Theorem I4.1U1 the root A = A(e) of 
E(X,e) = is infinitely smooth in e 1 / 2 . Therefore, the eigenvector (u,w) is also infinitely smooth in 
e 1//2 . By Lemma 14.11 the kernel of operators £+ and £_ is one-dimensional in L 2 (M), such that the 
limiting relation (|4,13|) follows by smoothness (eigenvectors are defined up to an arbitrary multiplicative 
factor) . 

(iii) It follows from the decay (|4.5|) and (|4.fi|) for the four fundamental solutions that there exist 
constants ^4±(A) and B±(X), such that 

M -> A± ( _^ \ e^ K - x + B ± (^ _ K + J e^+ x as x -> ±oo, (4.16) 

where k± are defined by the characteristic equation ()4.3j) with A = A(e). By expansion (|4.4|) . the 
expansion of slowly decaying solutions e^ K ~ x gives Q4.14JI . while expansion of fast decaying solutions 
e Tn+x [ s no t written. □ 
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Theorem 4.11 Let f(q) be C 2 (R+) and V(x) satisfy E3J). Let P^io + and S'\ vl0 / for the 
black soliton 4>o(x) of Definition \2. b\ Let A = A(e) be a small eigenvalues of the spectral problem \3.iy\) 
for sufficiently small e with ReA > 0. Then, the value of A(e) is given by the root of the characteristic 
equation: 

ReA > : (P' r \ vl0 ) A 2 - M '[ {So) A + eM"(s ) = 0(e 2 ). (4.17) 



r\vl0) 



Proof. By Remark 13, 161 the spectral problem p.lijj) for sufficiently small e can be written in the form: 

[L+ + eM + + 0(e 2 )] u = -Xw, [L_ + eM_ + 0(e 2 )] w = Xu. 
By Theorem 14.101 the eigenvector (u, w) and the eigenvalue A = A(e) can be expanded in powers of 

e l/2. 

u = <p' {x) + e 1/2 m + eu 2 + e 3/2 n 3 + 0(e 2 ), w = e l/2 w l + ew 2 + e 3/2 w 3 + 0(e 2 ), 
and A = e x / 2 Ai + eX 2 + e 3 ^ 2 A3 + 0(e 2 ). The first-order corrections terms satisfy the system 

L +Ul = 0, L-Wi = Ai$,. (4.18) 

By the expansion ()4.14|) . we are looking for a solution with u\ E L 2 (R) and w\ E L 00 (1R). By Lemma 
14. 1( the explicit solution is 

ui = ci0' o (x), wi = ai(j) (x) - XiImd v lI(x)\ v io, 

where (ai,ci) are parameters. Without loss of generality, we can set c\ = 0. The second-order 
corrections terms satisfy the system 

L+u 2 = -XiWi - M + (p' , L-w 2 = A 2 0' o (4.19) 

By the Fredholm alternative for L + , we obtain the constraint 

-X 1 (<t>' ,w 1 )-{<f>' ,M +( l> , ) = 0. 

By using (|3.14| and (|3.25|) . the constraint is equivalent to the characteristic equation 

i(P;Uo)A 2 + ^M"( So ) = 0. (4.20) 

Since w\ E L°°(1R) and the Fredholm constraint is satisfied, there exists a solution u 2 E -L°°(]R) of 
the first inhomogeneous equation ()4.19|) . There exists also a solution w 2 E L°°(M) of the second 
inhomogeneous equation (|4.19f) similarly to the solution w\ E -L°°(]R). By the expansion (|4,14|) . we 
shall add a homogeneous linearly growing solution (iii) of Lemma 14. II which is not in L°°(M), such that 
w 2 ^ L°°(R). As a result, the second-order corrections terms are written in the form 

u 2 = c 2 (j)' Q (x) + u 2 (x), 

w 2 = a 2 (f) (x) - X 2 Imd v U(x)\ vl0 + b 2 (x(p Q (x) -Imd v U(x)\ vl0 ), 

where (02,62,02) are parameters (c 2 = without loss of generality) and u 2 E L°°(R) is a solution of 
the inhomogeneous equation 

L + u 2 = -Aiai^o - M } S °\ nid v U(x)\ vi o - M + 4>' . 
K\vio 
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The third-order corrections terms satisfy the system 

L+U3 = -X2W1 - X\w 2 , L-W3 = A3 0o + Aiu 2 - M-W\. (4-21) 
By the Fredholm alternative for L+, we obtain the constraint 

-A 2 (0 o ,wi) - X 1 ((j)' ,w 2 ) = 
which is equivalent by virtue of (|4.2|) to the characteristic equation 

(^ko)AiA 2 + toW^ Xtb* = 0. (4.22) 

Combining two corrections Ai and A 2 in A(e) = e 1 / 2 Ai + eA 2 + 0(e 3 / 2 ), we rewrite the characteristic 
equations (|4.2U|) and (|4.22|) in the form 

(P;| 1 , i0 )A 2 + eq (S'\ vl0 )b 2 A + eM"(s ) = 0(e 2 ). (4.23) 
In order to find b 2 for ReA > 0, we need to consider w(x) = e l l 2 wi + ew 2 + 0(e 3 / 2 ) for large ±x S> 1: 
w(x) -> ±^o~ (e 1 / 2 ( Ql - AAG^o) +e{a 2 - X 2 8 v e ± \ vl0 + b 2 (x - ^O ± U i0 )) + 0(e 3 / 2 )) (4.24) 
where = lim 0(x) and lim 0o(x) = i^/go have been used. Matching the asymptotic expansions 

x— >±oo x— >±oo 

(|4.14() and (|4.24l) . we find a linear system on parameters (01, b 2 ): 

cb 2 = =fAi (ai - AAe^o) • 
The linear system has the explicit solution 

2cb 2 = xjd v (e+ - e-) \ vl0 , 2Aiai = a 2 ^ (e+ + e~) \ vl0 , 

such that 

, _ (S'\ vl0 ) 2 (S'\ vl0 )M"(s ) 
2c Al " 2c(PXio) ' 
where S 1 = + — 0~. As a result, the characteristic equation Q4.23[) reduces to the form (|4.17[) . □ 

Corollary 4.12 Let P/|„;o / and S"|^ /or t/ie fe'nfc mode of Theorems HOI fO^l 1^31 and 
\4-ll\ for sufficiently small e. Then, 

• If P' r \ v io > 0, the kink mode with M"(so) > /ias precisely one quartet of small complex eigen- 
values (N r = Nf = 0, N c = 1), while the kink mode with M"(sq) < has precisely one pair of 
small real eigenvalues (N r = 1, iV7~ = N c = 0) in the spectral problem \3.iy\) . 

• If P(.\ v iq < 0, the kink mode with M"(sq) > has precisely one pair of small real eigenvalues and 
one pair of finite real eigenvalues (N r = 2, = N c = 0), while the kink mode with M"{sq) < 
has precisely one pair of finite real eigenvalues and no small eigenvalues (N r = 1, Nf = N c = 0) 
in the spectral problem 13.191) . 
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Remark 4.13 Comparison of the characteristic equation (|4.17|) and the linearized version of the New- 
ton's particle law (|1.7|) shows that 



u ~P'\ A - 9o(y| {1 

A*0 — " r |i;J.O) A 



Ho J 



Both constants are positive if P/kiO > 0, i.e. if the kink is stable in the spectral problem (|3.1|) . 



Remark 4.14 Characteristic equation (|4.17j) can be derived from the power expansion of the Evans 
function E(X, e) of Definition 14.51 

E(X, e) = A (aA 2 + dA 3 + (3e + (5Xe + 0(A 4 , A 2 e, e 2 )) , (4.25) 

where (a,/3) are constants from the expansion (|4.15j) and (a,/3) are new constants. Computations of 
these constants from derivatives of E(X, e) are technically involved (see |19l I44| ). These computations 
are replaced in Theorem 14.111 with direct expansions of eigenvectors and eigenvalues of the spectral 
problem (|3.19j) in powers of e 1 / 2 . 



Example 4.15 Continuing Example 13.71 we consider the cubic NLS with f(s) = s and qo = 1, where 
P' r \vlo = 4 and S'l^o = 2. As a result, the characteristic equation (|4.17f) is written explicitly by 

ReA > : A 2 + -^M" (so) (l - - ) = 0(e 2 ). (4.26) 

This equation has only one real-valued root A(e) > for M"{sq) < and two complex-conjugate roots 
with ReA(e) > for M"(sq) > provided that e > is sufficiently small. The validity of the expansion 
()4.26|) will be tested in Section [5J 



Remark 4.16 If the characteristic equation (|4.26j) is formally applied to the cubic GP equation (|1.2j) 
with f(s) = s, qo = 1 and V{x) = x 2 , we obtain M"{s) = 2 J R sech. 2 (x)dx = 4, such that the 
characteristic equation (|4.26j) is 

A 2 -|A + e = 0(e 2 ). 

This characteristic equation was derived in [31] with a formal method for slow dynamics of dark solitons 
in parabolic potentials subject to radiative boundary conditions. The validity of radiative boundary 
conditions for parabolic potentials V(x) can not be verified by the present analysis. 



Remark 4.17 The characteristic equation (|4.17f) can be rewritten in the form 

(P/Uo) A 2 + go( ^° )2 A 3 = -eM"(s ) + 0(e 2 ). 

The left-hand-side of this equation was derived in Eqs. (19)-(20) of 35 and Eqs. (2.37)-(2.38) of |3Hj in 
a more general context of dark solitons with v S (— c, c). The method of 35^ was based on asymptotic 
theory for slow dynamics of dark solitons, while the method of [HE] was based on slow decay conditions 
for eigenfunctions of the linearized problem. Here we have replaced these formal methods with rigorous 
proof in the framework of Theorems 14.101 and 14.111 
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5 Numerical approximations of eigenvalues 



We test here the predictions of the characteristic equation 1)4.26)) for the cubic GP equation 1)1,2)) with 
f(s) = s and the two potentials V\(x) and V2(x) in 1)1.4)) . in particular, we focus on examining the 
dependencies of small unstable eigenvalues of the spectral problem ()3.19j) versus parameters e and n. 
The numerical approximations of the kink mode (j) e (x) are obtained by means of fixed point iterations 
of the ODE ()2.12|) . The iterations are applied to a finite-difference discretization of the computational 
domain x £ [-L, L] (typically L = 10) on a grid of N nodes (typically N = 1600) with a spacing Ax 
(typically Ax = 0.2). Subsequently, the spectral problem 1)3, 19|) is discretized in a matrix eigenvalue 
problem that, in turn, is solved through standard numerical linear algebra routines. 

In the case of the potential Vi(x), as is considered in Examples 12.151 and 13.19) the positive-definite 
sign of L"(0) (and hence M"(0)) leads to a sole kink mode bifurcating from so = (and staying at 
s = by Remark l2.10|) . The kink mode is located at the minimum of the effective potential M(s) and 
is unstable due to a complex quartet of eigenvalues according to the characteristic equation 1)4.26)1 . 

The numerical results on Figures |2] and El fully confirm the above picture. Fig. |2] shows only one 
solution (f> £ (x) of the ODE 1)2.12)1 with the potential V\{x) and a unique quartet of complex eigenvalues 
A = A r + iXi in the spectral problem (j3.19)) . The left panel of Fig. El shows the real part of this quartet 
as a function of k for a given e = 0.2, while the right panel shows the relevant real part as a function of 
e for a given k = 1. The predictions of the characteristic equation 1)4.26)1 are shown by dashed-dotted 
lines, while the numerically obtained eigenvalues are shown by thick lines. 

The non-monotonic behavior of the real part of complex eigenvalues is produced by the truncation of 
the computational domain x G [— L, L] and subsequent discretization on a finite grid. This numerical 
phenomenon is explained in (see their Figure 2) as follows. The continuous spectrum of the spectral 
problem ()3.19)) becomes a finite spectral band along the imaginary axis near A = due to the truncation 
and the band is represented by isolated eigenvalues due to the discretization. The quartet of eigenvalues 
bifurcates from the point A = in the direction of the imaginary axis with small real part for small 
e and interferes with eigenvalues from the discretized continuous spectrum. This interference leads to 
the non-monotonic behavior of the real part of the relevant eigenvalues on Fig. |21 We note that this 
effect is not present for real eigenvalues bifurcating from the point A = 0. 

In the case of the potential Vvix), we have a more interesting phenomenology. While the potential 
always has two maxima at x = ±2/k, the effective potential M(s) possesses two maxima at sq = is*, 
s* and a minimum at so = for k < kq ~ 3.21 and only one maximum at so = for k > hq. 
Branches of solutions with M"(sq) < are unstable due to a pair of real eigenvalues, while the branch 
of solutions with M"{sq) > is unstable due to a quartet of complex eigenvalues, according to the 
characteristic equation 1)4.26)1 . The transition of the kink mode s$ = from k < kq (when M"(0) > 0) 
to k > k,q (when M"(0) < 0) resembles the subcritical pitchfork bifurcation at k = k,q. 

The numerical results for the potential V^(x) are shown on Figures El El and H3 Branches of the 
solutions with sq = ±s* are denoted by a dashed line, the branch of the solution sq = with M"(sq) > 
is denoted by thick solid line and the same branch with M"(sq) < is denoted by thick dashed line. 
The corresponding predictions from the extremal points of the effective potential M(s) and from the 
characteristic equation l)4.26|) are shown by dash-dotted lines (thick for so = and thin for so = ±s*)- 

The bifurcation point is found numerically to be kq ~ 3.26 for e = 0.2 in a very good agreement with 
the value kq = 3.21 obtained from M(s) with M"(0) = 0. The computational error is approximately 
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Figure 2: Bifurcation results for the potential V\{x). The top panels show the kink mode (solid 
line) and the potential V\{x) amplified by a factor of 5 (dash-dotted line) for (e, k) = (0.2,1.1) (left) 
and (e, k) = (0.2,6.4) (right). The bottom panels show the corresponding spectrum of the linearized 
problem (the numerical result is shown by circles, while the theoretical prediction is shown by stars). 
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Figure 3: The left panel shows the theoretical (dash-dotted line) and numerical (solid line) dependence 
on k of the real part of the unstable complex eigenvalue for fixed e = 0.2. The right panel shows the 
dependence of the same quantity on e for fixed k = 1. 
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Figure 4: The subcritical pitchfork bifurcation in parameter k for the potential V2{x) for fixed e = 0.2. 
The left panel shows the center of mass sq of the kink modes (sq 7^ by dashed line, sq = by thick 
solid and dashed lines). The theoretical predictions of sq are shown by dash-dotted line. The vertical 
line gives the theoretical prediction for the bifurcation point k = kq. The right panel shows the real 
part of the unstable eigenvalues for the relevant kink modes, using the same symbolism as the left panel. 
The theoretical predictions of eigenvalues are shown by thick and thin dash-dotted lines, respectively 
for the branches with sq = and sq 7^ 0. 



1.5%. The values of so for the kink mode 4> t (x) are obtained numerically from its "center of mass" 
defined by 

f-^z(l-M 2 )<'x .... 

s » = JZJ^MW- (5 ' 1) 

The values of sq are plotted on the left panel of Fig. 01 for e = 0.2 in a good agreement with the value 
s* obtained from M(s) with M'(s*) = 0. 

The solution profile <^> e (x) and the corresponding linearization spectra for the different branches and 
for particular choices of (e, n) are shown on Fig. |SJ In agreement with the characteristic equation 
(|4.2fij) . the kink modes with so = ±s* for k < kq and with sq = for k > kq has a pair of small real 
eigenvalues, while the kink mode with sq = for k < kq has a quartet of small complex eigenvalues. 
The quartet of complex eigenvalues for the kink mode with so = and the pair of real eigenvalues for 
the kink mode with so 7^ exists for k < kq, merge at the origin at n = kq and lead to a pair of real 
eigenvalues for the kink mode with so = and k > kq . The right panel of Fig. 0] shows the real parts 
of unstable eigenvalues for each kink mode. 

Fig. H3 shows the dependence of the relevant eigenvalues versus e for a fixed k = 1 < Kq. In this case, 
three branches of kink modes exist and the branch with sq 7^ has a pair of real eigenvalues, while the 
branch with so = has a quartet of complex eigenvalues. We can see that the non-monotonic behavior 
of the real part of unstable eigenvalues is only observed for the quartet of complex eigenvalues. We can 
also see from all figures of this section that the agreement between numerical and theoretical results is 
excellent for e < 0.3 and deteriorate for e > 0.3 due to the truncation error of the order of 0(e 2 ) in the 
characteristic equation (|4.2fi|) . 
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Figure 5: The left quartet of panels shows the solutions with so 7^ (dashed lines) and the potential 
V2(x) (dash-dotted line) for (e, k) = (0.2,1.0) and (e, k) = (0.2,3.1). The corresponding spectrum 
features a pair of real eigenvalues (numerical results are shown by circles and the theoretical predictions 
are shown by stars). The right quartet of panels shows the similar picture for the kink mode with so = 
for (e, k) = (0.2, 1.025) and (e, k) = (0.2,3.975). The corresponding spectrum features either a quartet 
of complex eigenvalues or a pair of real eigenvalues. 
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Figure 6: The real part of the unstable eigenvalues versus e for k = 1 for the solution branches with 
sq ^ (left panel) and with sq = (right panel). The numerical results are shown by dashed and thick 
solid lines, while the theoretical predictions are shown by dash-dotted lines. 
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6 Numerical simulations of the GP equation 



We examine here the dynamics of the unstable kink modes in the full GP equation (|1.2|) by using direct 
numerical simulations. The time-evolution problem is approximated by the fourth-order Runge-Kutta 
method applied to the spatial discretization of the GP equation. The output of the fixed point iteration 
was used as input in the time evolution integrator with the time step At (typically At = 0.001). The 
results of the time evolution are compared against the effective Newton's particle equation (|1.7|) for 
the position s(t) of the center of dark soliton (f> e (x — s(t)). 

In order to test the theoretical result, we have to use the following numerical technique. The initial 
condition u(x, 0) of the GP equation (|1.2|) is specified in the form of the kink mode (j> e (x) plus a small 
(typically 10~ 4 ) perturbation multiple of its most unstable eigenmode. The time-evolution problem is 
integrated for an initial period < t < to, during which the dark soliton acquires a small speed due to 
instabilities, which quickly grows for t > to- At the time instance t = to, we approximate the values of 
s(to) and s(to) by using the center of mass (|5.1j) at t = to and at earlier time instances. The Newton's 
particle equation (|1.7j) is initialized at t = to with given s(to) and s(to) and then integrated with the 
fourth-order Runge-Kutta method. 

Figures[2|and|Hlillustrate the time evolution of an unstable dark soliton, which possesses a pair of real 
eigenvalues in the linearization spectrum. Fig. [7] corresponds to the potential V^fa) with k = 4 > kq, 
when the kink mode with sq = has a real eigenvalue A = 0.0241 (the theoretical prediction of the 
linearized Newton's particle equation is A ~ 0.0253). We observe from the figure that the unstable kink 
mode undertakes a monotonic transition to a stable dark soliton, which escapes the double-humped 
potential Vzfa) and travels with an asymptotically constant speed. The predictions of the Newton's 
particle equation shown by thick dash-dotted line captures the entire process accurately but slightly 
precedes the full time-evolution of the GP equation. The latter discrepancy can be attributed to the 
larger values of A for the unstable eigenvalues. 

Fig. |H] shows a similar monotonic transition for the potential V^(x) with k = 1, when the kink mode 
with so = s* ~ 2.23 has a pair of real unstable eigenvalues. In this case, the theoretical prediction 
A ~ 0.1251 exceeds the numerically obtained value A = 0.1211 too and the prediction of the Newton's 
particle equation precedes its counterpart from the GP equation. It is worth to note the qualitative 
agreement between the two time evolutions, including the small "leg" formed in the trajectory as the 
dark soliton passes the unstable kink mode with sq = — s* ~ —2.23. 

Figures l9l and HOl illustrate the time evolution of an unstable dark soliton, which possesses a quartet 
of complex eigenvalues in the linearization spectrum. Fig. |U] corresponds to the potential V^(x) with 
k = 1, when the kink mode with sq = has a quartet of complex eigenvalues with the real part Re A = 
0.06606 (the theoretical prediction of the linearized Newton's particle equation is ReA ~ 0.07724). We 
observe from the figure that the unstable kink mode oscillates in a local potential well of the double- 
humped potential V^(a;) with an increasing amplitude due to unstable complex eigenvalues. When the 
oscillations reach a large amplitude, the dark soliton escapes the maximum of the effective potential and 
transforms to a steadily moving soliton. The predictions of the Newton's particle equation represent 
this dynamics correctly with a larger deviation from the full GP equation in comparison with the case 
of monotonic transitions. 

Fig. 1101 shows a similar oscillatory behavior of an unstable kink mode with sq = for the potential 
V\(x) with k = 6.4. In this case, the theoretical prediction ReA ~ 0.0123 exceeds again the numeri- 
cally obtained value ReA = 0.0118 and the prediction of the Newton's particle equation precedes its 
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Figure 7: The unstable evolution of the kink mode with sq = for the potential V^(x) with k = 4 
and e = 0.2. The dashed-dotted line shows the result of the Newton's particle law initialized around 
t = 180. 
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Figure 8: The unstable evolution of the kink mode with sq ~ 2.23 for the potential V2(x) with k = 1 
and e = 0.2. The dashed-dotted line shows the result of the Newton's particle law initialized around 
t = 140. 
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Figure 9: The unstable evolution of the kink mode with sq = for the potential V^(x) with k = 1 and 
e = 0.2. The dashed-dotted line shows the result of the Newton's particle law initialized at t as 1000. 
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Figure 10: The unstable evolution of the kink mode with so = for the potential V\{x) with k = 6.4 
and e = 0.2. The dash-dotted line shows the result of the Newton's particle law initialized around 
t w 600. 
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counterpart from the GP equation. 

In the end, we mention that the rigorous derivation of the Newton's particle equation for slow 
dynamics of a bright soliton in an external potential has been reported recently in [6l ll3llTH] . Derivation 
of its counterpart (|1.7|) for slow dynamics of a dark soliton is an open problem of analysis. Our 
numerical results suggest that this Newton's particle equation is highly appropriate for understanding 
the nonlinear time-evolution of dark solitons in the GP equation with small external potentials. 

7 Conclusion 

We have systematically analyzed the persistence and stability of dark solitons in the presence of small 
decaying potentials. We have shown how the effective potential can be used to predict bifurcations 
of kink modes in a small potential and to approximate small unstable eigenvalues of the linearization 
spectrum. These theoretical results have been tested against the numerical bifurcation results indicating 
excellent qualitative and good quantitative agreement. We have also conjectured a dynamical evolution 
equation (the Newton's particle law) that can be used, quite successfully, to describe the motion of the 
kink modes and the manifestation of their instabilities. 

One of the directions of interest for future studies is to expand the present results to other types of 
potentials which include periodic and confining potentials. While, as argued in the text, we expect many 
of the qualitative features to persist, periodic or growing potentials may possess additional interesting 
properties due to the presence of spectral bands or purely discrete spectrum in the linearization problem. 

Another open direction would involve extending the present analysis to the two-dimensional setting 
and, in particular, to the case of vortices in the presence of external potentials. While some of the 
techniques applied herein (in particular, ones involving perturbative expansions) would apply to the 
latter case as well, others are more geared towards the one-dimensional setting (e.g. the Evans function 
technique). It would be especially interesting to generalize our current results to the two-dimensional 
GP equation. 
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